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THE  GENERATION  OF  AN  ACCURATE  AND  UNIFORM  TIME  SCALE 
WITH  CALIBRATIONS  AND  PREDICTION* 

Kazuyuki  Yoshimurat 

Atomic  Frequency  and  Time  Standards  Section 

National  Bureau  of  Standards 

Boulder,    Colorado     80302     USA 

Abstract 


We  express  a  predicted  time  interval  (or  frequency)  of  a 
single  clock  as  a  weighted  sum  of  frequency  data  obtained  by 
calibrations  against  a  primary  standard,    and  derive  a  matrix 
equation  for  the  optimum  weighting  coefficients   (called  the 
optimum  filter  response  function)  involving  the  Allan  variances, 
Two  approaches  are  used.     One  of  the  approaches  turns  out  to 
be  a  generalization  of  Barnes'  approach  described  in  his   1966 
IEEE  paper. 

We  solve  the  matrix  equation  to  get  the  optimum  filter 
response  functions  for  white  noise  frequency  modulation  (FM), 
flicker  noise  FM  and  linear  combinations  of  them.      Other 
important  time  dispersive  mechanisms  exist  in  practice  but 
are  not  considered  here.      We  obtain  the  result  that  the  mean 
square  time  prediction  error  -would  increase  as  elapsed  time 
t      for  the  case  without  intermediate  calibrations. 

We  obtain  the   result  that  with  a  small  number  of  good 
clocks  one  can  construct  a  time   scale  whose  accuracy  is 
limited  by  the  accuracy  of  a  primary  standard.     We  show  that, 
over  a  long  time  range,    linear  prediction  algorithms  based  on 
frequency  calibrations  with  a  primary  standard  give  a  time 
scale  of  much  better  accuracy  and  stability  than  when  inter- 
mittent calibrations  are  not  used,    and  that  (at  least  for 
statistically  identical  clocks),    no  advantage  is  gained  by  using 
the  time   scale  itself  as  a  "primary  standard"  for  intermediate 

calibrations. 

i 

Key  Words:  Accurate  and  uniform  time  scale;  Allan  variance; 
Dispersion  of  time  scale;  Ensemble  time  (error); 
Prediction  interval;  Primary  standard  and  clocks. 
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1.     INTRODUCTION 

An  atomic  time  scale  can  provide  one  of  the  best  time  scales  in 
terms  of  uniformity,    accuracy,    reliability,    and  reasonable  accessi- 
bility   [1].     In  the  areas  of  accuracy  and  uniformity,    it  has  a  clear 
advantage  over  astronomical  time  scales.     An  atomic  time  scale  can  be 
actually  realized  with  atomic  frequency  standards  and  appropriate 
counters;  such  devices  exhibit  undesirable  phenomena  such  as  frequency 
or  time  jumps,    aging,    and  failures,    as  well  as  the  signals  being  con- 
taminated by  random  noises.     All  of  those  phenomena  introduce  very 
important  problems  in  the  generation  of  an  atomic  time  scale.     Random 
noises,    however,    which  contaminate  the  signals  from  the  atomic  clocks-- 
commonly  white  noise  frequency  modulation  (FM)  and  flicker  noise  FM 
in  cesium  clocks--are  included  in  the  process  of  generating  a  time 
scale  and  result  in  time  dispersion  such  as  random  -walk  and  flicker  walk. 
So  these  random  noise  processes  can  be  said  to  contribute  basically  to 
the  time  scale  error,    and  are  the  only  processes  treated  in  this  paper. 

In  order  to  obtain  a  much  better  time  scale  than  that  generated  by 
a  single  clock,    a  group  of  clocks  is  commonly  used.     There  are  various 
algorithms  used,    such  as  a  weighted  sum  over  clocks  or  the  acceptance 
or  rejection  of  clocks,    to  get  as  good  a  time  scale  as  possible  from  the 
group  of  clocks.     These  algorithms  involve  assumptions  about  the 
behavior  of  each  clock.     In  practice,    laboratories  adopt  their  own 
algorithms  to  generate  local  time  scales  and  compare  these  time  scales 
with  each  other    [2]- [5].      It  is  common,    however,    for  laboratories  to 
try  to  obtain  a  more  uniform  time   scale  instead  of  a  more  accurate  and 
uniform  time   scale,    basically  by  taking  a  weighted  sum  (usually  0  or   1 
weights   [6])  over  several  clocks  without  using  a  primary  standard. 
This  can  result  in  a  significant  drift  away  from  the  defined  atomic  fre- 
quency and  an  accumulating  time  error  relative  to  an  ideal  time  scale   [7], 


In  a  recent  paper   [3],    Mungall  compares  and  discusses 
differences  between  time  scales  which  do  or  do  not  use  a  primary- 
standard,    and  points  out  the  problems  in  the  latter  time  scale.     It  was 
also  suggested  in  a  recent  conference  [5]  that  a  study  be  conducted  on 
how  to  have  concurrently  both  an  accurate  and  uniform  IAT  scale  (the 
International  Atomic  Time  scale).     In  order  to  accomplish  this,    a 
primary  standard  must  be  tightly  related  to  the  generation  of  the  time 
scale,    namely  it  must  be  used  to  calibrate  periodically  secondary 
standards   (or  clocks)  of  lower  accuracy.     Periodic  calibration  is 
usually  necessary  because  the  primary  standard  can  be  operated  only 
for  brief  periods  and  because  frequent  re-evaluations  of  its  accuracy 
are  necessary. 

In  this  paper,    we  predict  the  time  change  for  each  clock  between 
calibration  intervals  by  filtering  calibrated  frequency  data  and  derive  a 
matrix  equation  for  the  optimum  filter  response  functions  which  make 
the  mean  square  predicted  time  error  a  minimum.     We  apply  these 
results  to  the  generation  of  a  time  scale  from  several  clocks.     We  find 
that  we  can  predict  the  correct  time  with  quite  a  small  error  using  a 
small  number  of  clocks   (even  one),    and  thus  generate  a  time  scale  whose 
accuracy  is  limited  by  the  accuracy  of  the  primary  standard- -the  best 
case. 

The  dispersion  of  the  time  scale  thus  calculated  is  in  a  statistical 

sense  essentially  proportional  to  the  square  root  of  the  elapsed  time, 

while  the  time  scale  without  intermediate  calibrations  may  disperse  pro- 

i 

portionally  to  the  elapsed  time  for  the  flicker  noise  FM  case  which  is 
usually  the  dominant  noise  process  after  some  time   (more  than  10  days, 
for  example)  has  elapsed.     So  in  a  short  time  range,   the  latter  time  scale, 
which  may  use  a  large  number  of  clocks,   may  be  better  in  uniformity 
than  the  former;  but  in  a  very  long  range,   the  former  may  become  better 
and  better  in  accuracy  and  uniformity. 


An  approach  in  which  the  time  of  a  secondary  standard  is  predicted 
with  a  primary  standard  was  presented  in  Barnes'    1966  IEEE  paper  [8]. 
One  of  two  approaches  described  in  this  paper  is  a  generalization  of 
Barnes'  approach  to  include  any  number  of  calibrations  at  arbitrary 
positions. 

We  also  discuss  the  application  of  the  linear  prediction  method  to 
the  case  without  a  primary  standard,    trying  to  obtain  better  uniformity 
in  a  time  scale  than  that  due  to  accuracy  of  the  primary  standard  for 
both  short  and  long  time  range  (where  accuracy  here  is  used  for  only  the 
random  uncor related  contributions  of  the  primary  standard),    and  show 
that  (at  least  for  statistically  identical  clocks),    no  advantage  is  gained 
by  using  the  time  scale  itself  as  a  "primary  standard"  for  intermediate 
calibrations. 


2.     ALGORITHMS  TO  OBTAIN  THE  OPTIMUM  FILTER 
RESPONSE  FUNCTIONS 

Figure   1(a)  shows  an  ideal  case  in  which  the  time  of  a  clock  is 
moving  away  in  a  linear  manner  from  the  horizontal  line- -an  ideal  or 
defined  time  scale.     If  this  slope  is  determined,    the  whole  movement  of 
the  clock  times  can  be  perfectly  predicted,    so  we  can  construct  an  ideal 
time  scale  using  this  clock  by  subtracting  a  linear  term  from  its  time. 
But  signals  from  actual  clocks  are  contaminated  by  noise- -commonly 

by  white  noise  FM  (spectral  density  of  frequency  fluctuation   S    (f)  cc  f   ) 

-  1  y 

or  flicker  noise  FM    (S    (f)  cc  f      )or  both  (referred  to  symbolically  as 

0-1  y 

f     and  f        FM  noise).     In  this  case  we  cannot  predict  perfectly  the 

movement  of  the  clock  time,   thus  leaving  some  error  of  prediction  as 
shown  in  figure   1(b).     But  by  making  the  predicted  error  as  small  as 
possible,    one  may  construct  a  time  scale  close  to  an  ideal  time  scale. 
By  taking  a  weighted  sum  of  prediction  errors  for  a  group  of  clocks,   we 
may  obtain  a  smaller  error  in  the  time  scale  as  shown  with  the  dotted 
curve  in  figure   1(c).     We  call  this  weighted  error  the  ensemble  time 
error,    or  (for  simplicity)  the  ensemble  time.      The  ensemble  time  thus 
calculated  can  be  used  as  a  time  reference  (taking  the  place  of  an  ideal 
time   scale)  with  respect  to  which  the  time  of  any  clock  can  be  deter- 
mined. 

We  will  describe  two  approaches,    or  two  prediction  methods, 
using  frequency  data  for  a  single  clock  obtained  by  calibrations  against 
a  perfect  primary  frequency  standard.      Deterministic  frequency  offset 
(or  time  drift)  and  frequency  drift,    however,    are*  not  considered  here. 
Frequency  or  time  jumps  and  any  kinds   of  clock  aging,    -which  are  also 
important  phenomena  to  be  considered  for  the  generation  of  a  time  scale, 
are  not  dealt  with  either.      With  no  loss  of  generality,    we  will  assume 
that  the  noise  processes  are   stationary  through  all  the  calculations.      We 
must  determine  how  best  to  predict  the  time   slope  or  frequency. 


A.      Approach  using  calibrated  frequency-  data  which  are  located 
before  a  prediction  interval 

The  first  approach  assumes  the  calibration  intervals  precede  the 
prediction  interval  as  shown  in  figure  2.     Let  the  duration  of  the  predic 
tion  interval  be    T    and  the  duration  of  the  calibration  interval  be   hT. 
Here,    hT  and  the    It    are  not  necessarily  integral  multiples  of    T. 


Over  a  frequency  calibration  with  a  primary  standard,    the 
ge  fractional  frequen 
only  to  noise  is  written  as 


th 
average  fractional  frequency  offset  on  the     i       calibration  interval  due 


t 


y.  =  rr  /  y(t)dt  = 


*/ 


x(t^)    -  x(t^  -   hT) 


I  hT     /    '  w  hT 

t    -hT  (1 

I 

=     (x    -  x„      ,  )/hT 


where   x   is  the  clock  time  referred  to  an  ideal  time  scale.     Let  us 
define  the  predicted  average  fractional  frequency  over  the  prediction 
interval  as  a  summation  of  eq  (1): 


z2w    Yl 


y0  =  z.a«v     2^a«= :  (2: 

Xj    —    Xj  _  Xj    —    Xj 


where  the    a     are  weighting  factors  for    y    .      That  is,    the  numbers    a 

Xj  Xj  Xj 

weight  the  various  measurements  made  in  time  on  the  single  clock  bein^ 
considered.      In  a  very  real  sense  the    a      define  a  response  function  in 
the  time  domain  for  a  prediction  filter. 


We  define  the  predicted  time   x      at  the  end  of  the  prediction 
interval  as    x     =   x     +  Ty    .      Then  the  error  is    e     -    x     -  x   .      Using  (2), 
we  can  obtain  then 


eo    =   TV<x0-xi)  =  ^?a*Ax*,h"  Axo,i  (3) 


where    Ax„  ,    =  x. -  x_      ,      and      Ax„    ,    =   x^   -  x    . 
je,  h  A        l+h  0,1  0  1 

The  mean  square  prediction  error  is  written  as  (4) 


«o>  =  \  [?a.ux(hT)  +  2SE,  v.'  {-ux(kT)  +  H«k+  h)T> 


+    |ux((k-h)r)]l    +   "^Ea     (U(iT)  +   u((A+h-l)T) 


h       I 


(4) 


U   (u+h)r)  -  U  ({I  -  1)t)}  +    u  (t) 

X  X  X 


2  2 

where     (c_)    denotes  the  ensemble  average  of    e    ,    and   k    =     #    -    I     and 

U  (kr)  =     2[R    (0)  -  R    (kr)l     (see  reference   [9]).     R    (T)    is  the  autocor- 

XXX  X 

relation  of  time. 

Using  the  relation  between  the  U-function  and  the  Allan  variance 
with  zero  dead  time    [9  ]  , 


U   (kT)   =     -k(k-  1)  T2<CT2(k,  T)>   +  k     U   (T),  (5) 

x  y  x 


we  can  obtain 


H  Jo  ^  JS  x»  ^ 

*  (6) 

where    C(k)    and   Q(#)    are  functions  of  the  Allan  variances  given  by 

C(k)   =    -k(k-  l)(a2(k,  T)>    +     j(k+  h)(k+h  -    1)  <(T2(k  +  h,  T)> 

y  ^  y 

+    |(k-  h)(k-  h  -   l)<cr2(k  -  h,  T)> 

2  y 

and  .  (7) 

Q(£)  =  4  C- *(i  -   1)  (uZ(A,  T)>  -  (£+h  -  l)(£+h-  2)  <CT2(£+h-  1,  t)> 
2  y  y 

2  2 

+     (j&+h)(j&+h-l)<Cy    (£+  h,  T)>  +  (1  -  1)(1  -  2)  <a    (I  -  1,  T)>}. 


Equation  (6)  is  the  expression  for  the  mean  square  of  the 
predicted  time  error  in  terms  of  the  Allan  variances,    frequency  stability 
measures  in  the  time  domain. 

In  order  to  obtain  the  minimum  value  of  eq  (6),    let  us  consider  a 
function 


where    \   is  a  Lagrange  undetermined  multiplier.     Differentiating    F   with 
respect  to    a      with     I   =    I      or     £   =    I     and  canceling  out   X.    we  get  a 
matrix  equation  for  determining  the  optimum  filter  response  functions; 


[C(4    -    ^)-H}a       +     {H-C(l    -    ij)}a 

1  p 


Ea    fC(£-    £  )  -C(£  -  £J,}    =  h{Q(£  )  -  Q(£.)} 
£  p  1    J  p  1 

£^£ 
P 

(£     =    £       £       .  ..,    £    )  (8) 

p  Z       3  L 


*l 


and    T*    a     =   1     where    H  s    h(h-  1)  <cr2(h,  T)>. 


£  y 

=  £, 


In  eqs   (6)  to  (8),    h  and  the    #   need  not  necessarily  be  integers  as 
we  will  see  later;  the  Allan  variances  in  such  a  case  may  be  interpolated. 
The  Allan  variances  for  the  various  noise  models  are  calculated  and 
given  in  reference  [10]. 

First,    consider  the  white  noise  FM  case  where  the  Allan  variance 
takes  a  simple  form  because  of  its  independence  of  k: 

<CT2(k;  T))    =    CT2(T),         C(k)    =  hV(T),        Q(£)    =  h<J2(T), 

y    V        >  //  y  y  y 

2  ,2 

where     0"    (t)    =     \0"    (2,  T) ) .      From  the  matrix  equation,    we  get  constant 

y  y 

factors  for  all    £,a      =   1/L,    so  the  normalized  mean  square  error  will 

AJ 

be: 


V  1V     2        ,  _1_ 

2    2,         =    h^a£  +    l     =     Lh 


T   ay(r)  £  (9) 


1 

1 

as 

Lh  »  1 

2 

as 

=   1 

1 

T    'U 

as 

<<  1    . 

Notice  that  the  optimum   a      does  not  include  the  clock  parameter 

Ju 

2 
(j    (t),    but  only  the  parameter    L,   the  number  of  calibrations  used  for 

prediction.      For  large  Lh,    eq  (9)  approaches   1,   the  smallest  value. 

For  flicker  noise  FM,    the  Allan  variance  takes  a  more 

complicated  form,    but   C(k)    or    Q(&)    can  be  written  in  a  product  form  of 

the  Allan  variance  and  other  parameters. 


(1)  k  -    1       ' 


so 


C(k)    =    a2(T)c(k)/c(l),        Q(j>)   =    CT2(T)q(i)/c(l) 

y  y 


where 


c(l)      =      2   In  2 

c(k)      =    -k2in(k)  +  ^(k  +  h)2jen(k  +  h)  +^(k  -  h)2j^n(k  -  h)  (10) 

1  2  2 

q(j&)      =^{4    ln(l)-(i+h-l)     £n(jfc+h  -    1) 

2  2 

+  (A  i-h)    ln(l+  h)  +  (&  -   1)    £n(JL  -   1)}  . 


The  mean  square  error  and  the  matrix  equation  will  be  given  by 
,2 

2  2°  -  =  -5p>*vc<k)+!?a.eq(fl-?a^n,h,      (u 

T  a    (T)/c(l)  h      £<  A'  4  I 


and 


[c(j>     -    j&)-H}a      +    {H  -  c(£     -    X.)}a 
p         1  j&.  p         1        £ 

r  1  r  p 


A      '   P 


I)  -  c(i-    A  )}==h{q(X   )   -  q(l.)] 

1  p  1 


4^4 


1 


(  x     -     Jb   ,     H   ,    . 
p  12 


(12) 


respectively.     Notice  that  the  matrix  equation  does  not  contain  the  Allan 

variance,    so  the    a     depend  only  on  the  calibration  parameters. 

0         -1 
For    f     +  f        FM  noise,   the  Allan  variance  is  given  by  the  sum  of 

the  Allan  variance  of  each  process,    assuming  the  processes  are 

independent; 


<CT    (k,  T)>      =       <<72(k,  T)>  +       <CT2(k,  T)>       . 

7  f  Y  f" 


So  the  mean  square  error  is  given  by 


<eo> 


2    2,    , 
T    (T    (T) 

y    f-i 


=    r 


1^2       , 


c(l 


"T  EEa    a    ,c(k) 
L     h2   *<*'     £   X 


(13) 


and  the  matrix  equation  by 


10 


[c(JL     -    I)  -H}aA     +    {H  -  c(i     -    I  )}a 
P         1  Aj  P         1 


P 


+   V*  a   {c(£-    A  )-c(£-    1J}  =   h{q(X  )-  q(l.)} 

*L-/      X  p  1  pi 

X  7^  J^ 

P  (1     =  £_,    JL,    ....    jL)  (14) 

p  1  Z  Li 

2  2  2 

where    r   =  r(T)   =   CT    (T)   „/cr    (T)      ,    and   H  =  h    4n(h)   -   rhc(l). 

y     fo     y     f-i 

Figure  3  is  a  computed  result  for  a  prediction  with  only  one 
calibration  interval.      The  horizontal  axis   shows  the  time  distance  of 
the  calibration  interval  from  the  prediction  interval.     The  vertical  axis 
shows  the  mean  square  prediction  error  normalized  by  that  with  cali- 
bration position  equal  to  1.      The  prediction  error  of  the  white  noise  FM 
case  does  not  depend  on  the  position  of  the  calibration,    while  the  pre- 
diction error  of  the  flicker  noise  FM  case  is  getting  worse  with  the 
distance.     So  it  is  better  with  flicker  noise  FM  to  locate  a  calibration 
interval  as  close  to  the  prediction  interval  as  possible. 

Figure  4  is  the  case  where  from  two  to  five  calibration  intervals 
are  used  for  a  prediction.      The  last  calibration  is  always  located  just 
preceding  the  prediction  interval.      But  positions  of  others  are  change- 
able,   which  are  shown  on  the  horizontal  axis  as  "spacing.  "     The  vertical 

2 
axis   shows  the  mean  square  prediction  error  normalized  by    T      times 

the  Allan  variance. 

As  expected,    with  white  noise  FM,    the  prediction  error  does  not 

depend  on  the  position  of  calibration,    only  on  how  many  calibrations  are 

used  for  a  prediction  (actually  on  the  total  calibration  time).      With 

flicker  noise  FM,    it  depends  on  the  position,    as  well  as  the  number  of 

calibrations  used.      In  the  two  cases  for  flicker  noise  FM,    there  are 

optimum  spacings  which  are   shown  in  the  diagram  with  black  points  where 

the  mean  square  error  becomes  minimum. 
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Figure  5   shows  the  optimum  filter  response  functions  with   L 
equal  to  5.     One  is  for  a  two -day  (t   =  1  day)  spacing  between  calibrations 
and  the  other  for  contiguous  calibrations.      For  white  noise  FM,     a 
factors  are  always  constant,    while  for  flicker  noise  FM,    the  calibration 
interval  nearest  the  prediction  interval  is  always  given  the  largest  value, 
showing  that  it  may  be  enough  to  use  only  one  calibration  interval  nearest 
to  the  prediction  interval.     With  contiguous  calibration,   the   a     are 
oscillating,    perhaps  due  to  an  aliasing  phenomenon. 

B.      Approach  using  calibrated  frequency  data  which  are  located 
in  a  prediction  interval 

In  this  approach  using  calibrated  frequency  data,    the  calibration 
intervals  are  located  within  the  prediction  interval  as  shown  in  figure  6. 
M    and   I  are  not  necessarily  intergers.      The  calculation  procedure  is 
almost  the  same  as  that  described  for  approach  A. 

The  average  fractional  frequency  offset  due  to  noise  over  the 
calibration  intervals  is  given  by 


H=  (Vx*-i)/T        <*  =  V  V  ••••  V-         (15: 


We  define  an  estimated  average  fractional  frequency  offset    y-   ,  , 

&  OM 

during  the    0  -  M    interval  by  a  summation  of  (15),    and  a  predicted  time 
interval    (x      -  x    )    as  follows:  , 


*L 


y0M 


Za„  y„  ,        V^  a„  =    !>      K,  -  x^   =    MTy\    _.  (16) 

I  y  I  L~i     I  MO  OM 


1  =  1 
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The  predicted  time  error  is 


€0M    =    "M"XM=    MT%  "   (XM"X0) 


=    Mj"  a„  Ax,    -   Ax„„ 
<-*     I        I  OM 


(17) 


where    Ax„    =    Ty„    =  x„   -  x„     ,     and    Ax„,,  =   x,  ,  -  x^. 
I  yi  l         JL-l  OM  M        0 

Then  the  mean  square  of  the  predicted  time  error  is  finally  given 

by 


<€        >  =    -2M    T2EE  a   a    ,C(k)  +  2Mt2£  a  Q(£)  -  T  M(M  -  1)  (a    (M,r)> 

(18) 
where   k    =    I     -I,    and 


C(k)  =    -k(k-  l)<CT2(k,  t)>  +    4(k+1)k<^(k+  1»  T)> 
y  2  y 

+    -j(k-  l)(k-2)  <a2(k-  1,  T)> 
2  y 

and  (18') 


Q(£)  =  4  {-(M  -  £)(M  -JL-l)  <a2(M  -  I,  r)  > 

2  y 


+     (M  -  4  +  1)(M  -  1)  (cr2(M  .  jj  +  I,  T)>   -   (X  -  1)(A  -  2)  <CT2  (£  -1,7)) 


+     Jt(i  -  1)<CT2(^  r)>] 

y 


13 


From  eq  (18),    the  matrix  equation  for  the  optimum   a     is 


C(l .-   1,).      -C(i    -   i1,ai 

1  P 


£>£{cu-  y  -  c(i-  x1)3  =  ^{Q(Ap)  -  Q(^)} 


i^i1 


(i  =  v  v  •••'  V  <19> 


and  V*  a  .  =    1. 

i  * 

For  white  noise  FM,   the  optimum  response  functions  again 
become  constant.     Then  the  normalized  mean  square  error  will  be  given 
by 

<   2      ) 

OM  .,2^     2       w  M  -  L 

-_-    =    M    E  a,   -  M     =    -j— 

7    CTy(T)  *  (20) 

2 

-  M    /L.  as      M  »L 

-  0  as      M    =    L 

-  M(M  -  1)    as      L    =     1 

where    a      =     1/L. 

For  flicker  noise  FM,    the  mean  square  error  is  given  by  eq  (21) 
and  the  matrix  equation  for  the  optimum   a      by  eq  (22). 

Xj 

,2     . 

=     -2M2X)Ea/,a/),  c(k)  +    2MX;a.q(£)   -  M2Xn(M)  .        (21) 


T     (T    (T)/c(l)  £<£  £ 

y 
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c(£     -    £  )a       -  c(£     -    £  )a 
p         1      £.  p         1      £ 

1  p 


+  Z)a^[cU_  y  ■ c(i_  v]  =  m  [q(y  -  q(v] 

£tl, 
£^£ 

P  (4     =  I       I       ....    I    )  (22) 

p         2        3  L 


where 


c(k)   =  -k2£n(k)  +  ^(k  +  l)2£n(k  +   1)  +  hk  -    l)2£n(k  -    1) 

and  (23) 

12  2 

q(jfc)   =    -  {-(M  -  I)     4n(M  -  £)  +  (M  -  £  +  1)     £n(M  -  4  +   1) 

-    (£   -    l)2  4n(4  -    1)  +     £2  £n(£)}     . 


Now  let  us  consider  the  case  where  there  is  only  one  calibration 
in  a  prediction  interval  in  figure  6.      Equation  (23)  will  be  reduced  to 


where 


<€2      ) 

—TT~     =    "T77   (2q(Jt.)   -  M  4n(M)}  (24) 

c    c.  c  (1  )  1 

T    CT    (T) 

y 


q(0   =    -r  {-(M  -  £   )2je.n(M  -  £  )  +   (M  -  1   +   l)2£n(M  -  £    +  1) 

x  Cj  ±  J.  A  J. 

-   (I    -    l)2£n(£    -    1)+    A    £n(£  )]    . 
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In  eq  (24),    for  example,    let   M  =  3    and  locate  the  calibration 
interval  in  the  middle  of  the  prediction  interval.      Then  the  result 
coincides  with  eq  (31)  of  Barnes'  paper  published  by  IEEE  in  1966  [8]. 
And  for  another  example,    let   M  =  2n  +  1,    where    n   is  an  integer,    and 
again  locate  the  calibration  interval  in  the  middle  of  the  prediction 
interval- -namely,     I     -  n  +  1;    the  result  again  coincides  with  eq  (61)  of 
Barnes'  paper.      Thus  eq  (24),    where    I     and   M   are  not  necessarily 
integers,    seems  to  be  the  analytical  solution  of  Barnes'  eq  (47).     Thus 
approach  B,    which  can  include  any  number  of  calibrations,    can  be  said 


to  be  the  generalization  of  Barnes'  approach. 

For    f     +  f         FM  noise,   the  meai 
equations  will  be  given  by  (25)  and  (2  6) 


For    f     +  f         FM  noise,   the  mean  square  error  and  matrix 


<4> 


2    2 

T    (J    (T) 

y    r1 


■^<£*rU)  +  W)  H2EEa^,c,k, 


l<l 


+    ZMj^a    q(l)   -  MZin(M)l 


(25; 


2  ,    2 

where    r   =  r(T)    =   a    (T)  _/cr    (T)       , 

y     fo     y     f-i 


[rc(l)  +  c(£    -  j^)}  a       -{rc(l)  +  c(j&    -    J^)^ 

.    1  P 


+  2a£fc(1-  t   >-c^-  V]   =M[qUp)_  q(V} 
14  I, 


a  4  il 
v 


(l     =    JL,    l.,    ...,    I    )  (26) 

p  2        3  L 
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Figure  7  shows  the  mean  square  prediction  error  for  one 
calibration  in  a  week.      Calibration  duration  is   supposedly  one  day. 
The  horizontal  axis  is  the  position  of  the  calibration  interval.      With 
white  noise  FM,    the  mean  square  error  does  not  depend  on  the  position. 
With  flicker  noise  FM,    locating  a  calibration  in  the  center  gives  the 
best  result.     It  can  be  shown  that  predictions  based  on  calibrations 
located  symmetrically  about  the  center  of  the  prediction  interval  give 
the  same  error. 

Figure  8  shows  the  mean  square  error  versus  number  of 
calibrations  in  a  week.      Again  calibration  duration  is   supposedly  one 
day.      For  flicker  noise  FM,    optimum  calibration  distributions  are 
selected  for  calibration  number   1  to  3.      Quick  improvement  with  the 
number  of  calibrations  is   shown.      The  dotted  curve   shows  the  improve- 
ment due  to  a  number  of  clocks  without  calibration,    assuming  each  has 
the  same  noise  level,    and  assuming  constant  weighting  factors  for  each 
clock. 

Figure  9   shows  the  influence  of  calibration  duration.      Two 
calibrations  per  week  and  optimum- -for   f       noise  FM- -locations  for 
calibrations  are   selected.      The  vertical  axis  is  the  mean  square  error 
normalized  by  that  of  the  24-hour  calibration  duration.      As  expected, 
longer  calibration  duration  is  preferable. 

Figure   10   shows  the  optimum  filter  response  functions    a      for 

Xj 

flicker  noise   FM.      It  is  interesting  to  note  that  two    a     factors  for 

Xj 

calibrations   symmetrically  located  have  the  same  value,    analogous  to 
the  white  noise  FM  case. 

Figure   11   shows  the  influence  of  increasing  the  size  of  a  single 
prediction  interval,    M    ,    where  a  rate  of  two  calibrations  per  week  is 
used  and  their  positions  are  fixed.      For  both  white  noise  FM  and  flicker 
noise  FM,    the  error  is  almost  linearly  increasing  with  the  size  of  the 
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2 

prediction  interval,    instead  of  increasing  with   M      for  flicker  noise  FM. 

(Similar  results  were  obtained  by  Barnes   [8].  )     The  accumulated  mean 
square  error  over  M/7  prediction  periods  of  seven  days  each  lies 
essentially  on  the  same  curve;  we  shall  discuss  this  in  detail  later. 


3.     ENSEMBLE  TIME  USING  SEVERAL  CLOCKS 

With  approach  B,    it  may  be  possible  to  increase  the  size  of  a 
prediction  interval  and  construct  a  useful  ensemble  time  with  only  one 


P 


rediction  interval.     From  the  predicted  time   x        of  each  clock  at   M, 


M 


we  construct  an  ensemble  time  estimate  using  weights    w^  ,    for  each  of 

b  to  M 

the   N   clocks: 


N        .        .  N 

-ens         <e^      i     „i         ^->      i 

M  X— i      M    M       L*i      M 

i=  1  i  =  1 

The  correct  value  for  this  ensemble  time  is 


N 
mens         v-^*      i       i 

T  =      >     w      x 

M  L~i     M    M 

i  =  1 


so  the  ensemble  time  error    I,    is 


or  \a 


*  Tens  *ens 

M  M        "       M 


N 

i       i 


^     =    I>M«M    "  (2?) 

i  =  1 


The  mean  square  of   T       in  eq   (27)  is 


<*£>  =  ?>£<&>  <28> 

l 

since  the  clocks  are  independent. 

The  optimum  value  of  the  weighting  factors  which  makes  the  mean 
square  of  the  ensemble  time  error  minimum  will  be 

OM7  0M; 

by  the   same  calculation  procedure  as  described  for  eq  (8).      The 
minimum  value  of  the  mean  square  error  is  then: 


(?*) 


(^5>     ■      : 
M  min 


OM 


For  white  noise  FM  or  flicker  noise  FM,    the  mean  square  error 
of  one  clock  involves  only  the  Allan  variance  to  describe  the  clock.      So 
under  the  assumption  that  all  the   N   clocks  have  the  same  type  of  noise- 
white  noise  FM  or  flicker  noise  FM,    the  optimum  weighting  factors  and 
the  minimum  mean  square  error  will  be  expressed  in  simple  forms  as 

opt      \7  a'W     al2(r) 

y  y 


and 

-1 


(FTfe) 


<*£>      ■      =    K.-T^y    — i—  (32) 

M  mm  M 


y 

where    £L       is  constant,    calculable  from  eqs   (20),    (21),    and  (30) 
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In  general,    for  mixed  noises,    optimum  weighting  factors  can  not 
be  expressed  only  with  the  Allan  variances,    but  include  other  parameters 
related  to  calibration  and  prediction  methods. 

For  the  general  case,    with  either  approach  A  or  B,    where  many- 
prediction  intervals  are  connected  in  series   (see  fig.    19),    the  change  in 
ensemble  time  over  the  (M  -    1,  M)  interval  can  be  written; 

N 
mens      mens  ^-*»       i  ^-s      i 

TM     "TM-1     =      2^WM(XM-XM-1)'LWM  =  1' 
i=  1 

and  its   estimate 

N 
-ens       -ens  <r-^      l     *  „ 

T  -T  =       >     w      (x       -  x  ). 

M  M-l  Z-r      Mv   M         M-l 

i  =  1 

Then  the  change  in  ensemble  time  error  is 

N 


i     i 

M       "M-l         £-*  "MCM     ' 
i  =  1 


!    i =  Z v ; :' 


or 

M     N 


j=li=l 


(TQ   =    0)   ,  (33) 


where  the    c.    is  the  prediction  error  for  the  i       clock  on  the   (j  -  1,  j) 

J 
interval,    and  the  initial  error  is  assumed  to  be  zero. 
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The  mean  square  of  eq  (33)  will  be 


<f2> 
M 


i  J 


EiZ,iZ,  v-v^     1     x     /    x    1  \ 

W.       {€.      >      +        22-r2-rW-W-'^e-'C-'/ 


i    Lj        J        J 


J<J 


J     J         J    J 


(34) 


under  independence  of  each  clock. 

The  mean  square  ensemble  error  now  contains  covariance  terms, 
in  addition  to  the  mean  square  terms  we  have  had  before.      Equations 
for  the  covariance  terms  can  also  be  derived  by  calculation  procedure 
similar  to  those  already  used. 

We  can  also  construct  a  matrix  equation  for  weighting  factors 
from  eq  (34) 


i2 


2w.  <€.    )  +    2  Y]  w.,    <€.-€.,)  +  \.   =  0 

j     J  frf.     J        J     J  J 

j  H 

(i  =  1,    2,    .  .  .  ,    N    and  j  =  1,2,    .  .  .  ,  M] 


(35; 


and    2J  w.   =  1    for  all   j    where    \.    is  a  Lagrange  undetermined  multiplier, 
i       J  J 

But  the   size  of  this  matrix  is  the   square  of    N   times    M,    which  will 
unlimitedly  increase  with  M,    and  the  -weighting  factors  will  be  functions 
of  M.     This  seems  impractical. 

For  white  noise  FM,    however,    since  the  covariance  term 
disappears,    the  results  become   simple. 


W.opt         U^     ,12,  ,   12,         (?      12.    J  12.    . 

J  \  1       <€.     >/         <€.     >  V  l      a       (T)/  CT      (T) 


(36) 
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and 


<*M>  = 

min 


In  general,    certain  clocks  in  the  ensemble  of  clocks  may  come  to 
be  rejected  with  increasing  time  on  account  of  their  prediction  errors 
getting  much  larger  than  others.     Then  it  may  be  preferable  to  construct 
the  ensemble  of  clocks  with  a  smaller  number  of  clocks  all  of  which  have 
the  same  type  of  noise.     Approach  B,    in  this  sense,    seems  better  than 
approach  A  because  it  may  give  a  simple  form  for  the  weighting  factors 
even  for  mixed  noises   (see  fig.    11). 

Under  the  assumption  of  a  fixed  weighting  factor  for  each  interval 
(w.   =  w  ),    again  we  can  get  simple  forms: 


<T^>     =£wi2<(lV)2>       .  (38) 

i  J 


So 


-1 
1  \  1 


i     ((Ec1)2)/      ((Ec1)2) 
J    J  j    J 


(39) 


and 


(40) 


We  calculated     ((Ic  )    )    in  (38)  in  order  to  obtain  some  numerical  idea  of 
the  mean  square  error  of  the  ensemble  time.     For  approach  B,    we  already 
showed  an  example  of  this  quantity  in  figure   11,    so  here  some  examples 
for  approach  A  will  be  explained. 
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Figures   12  and  13  show  the  influence  of  positions  of  calibrations. 
An  almost  equally  spaced  distribution  will  give  the  best  result  for  the 
flicker  case  as   shown  with  (1).      Figure   12  is  for  one  calibration  used 
for  a  prediction  and  figure   13  is  for  two  calibrations  used  for  a  predic- 
tion.    It  should  be  pointed  out  that  the   slope  of  the  mean  square  errors 

*  2 

T    ,(M)  may  actually  be   smaller  for  flicker  noise  than  for  white  noise, 

although  this  is  not  indicated  in  the  figures,    where  the  normalization 

influences  the  slope. 

Figure   14  is  the  accumulated  error  in  one  week  versus  the 
number  of  calibrations  in  a  week.     Rapid  improvement  ■with  increasing 
rate  is   shown  for  flicker  noise  FM. 

Figure   15  is  the  accumulated  error  versus  the  number  of 
calibrations  used  for  a  prediction.      Calibration  rate  is  two  calibrations 
per  week.      For  -white  noise  FM,    rapid  improvement  with   L    is   shown 
and  for  flicker  noise  FM,     L   =   1    gives  a  local  minimum  independent  of 
M.      For    f    +  f        FM  noise  with    r   =   10,    it  is  apparent  that  there  are 
optimum  values    L.     Improvement  with  larger  L,    however,    seems  quite 
small,    so   setting  L  equal  to   1  may  be   sufficient  for  flicker  noise  FM  or 
white  and  flicker  FM  noise  in  approach  A. 

Figure   18  shows  the  simulation  results  of  a  time   scale  using  one 
clock  for  mixed  noises   shown  in  figures    16  and   17,    using  a  white  noise 
and  flicker  noise  model   [11].      In  the  diagram,    improvement  by  pre- 
diction is  clearly  shown  with  curves  A  and  B,    which  are  predicted  time 
errors  for  curve  S,    the  time  error  of  the  mixed  frequency  noise  given 
in  figures   16  and   17.      We  will  discuss  figure   18  in  more  detail  later. 
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4.     CONSIDERATION  OF  "ACCURACY" 

All  the  calculations  so  far  mentioned  took  into  consideration  only 
noises  of  clocks- -white  or  flicker  noise  FM.     But  actually,   precision  of 
measurements  and  "accuracy"  of  the  primary  standard  must  be  also 
considered.     Here,    effects  of  random  uncorrelated  noise  as  it  affects 
the  accuracy  of  the  primary  standard,    denoted  "accuracy,  "  will  be  dis- 
cussed;  the  imprecision  of  the  measurements  will  be  assumed  negligible, 

We  shall  assume  that  the  total  average  fractional  frequency  offset 
may  be  expressed  as  a  sum  of  two  independent  random  variables, 


— i  — Oi  —  , ..  , . 

'a  =  yt  +  6yt.  (41) 

J  J  J 


where  the  first  term  depends  on  clock  noise  and  the  second  on  noise 
(inaccuracy)  in  the  primary  standard.     Then  the  predicted  average 
fractional  frequency  offset  will  be 


Ai        --,    i  _i  _    i  — Oi        _  ■  i       —  .    . 


and  the  predicted  time  error  will  be 


i  ^i         A    i 

«j     =    Tyi    -   Axj  (43) 


Oi    .^     i   .       -    .  Oi   .      _    i 

J  J 


=   «?i+?V(T«V^"+   ^ 


J 


where 


<j  -  EVTyi  >  -  Axj  and  6ej  s  £ai.(T6yi.' 


(44) 
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Then  the  mean  square  prediction  error,  finally,  must  be  written 
as  a  sum  of  the  noise  dependent  term  and  an  "accuracy"  dependent  one 
as 

,    i2.  ,    0i2.  ,  „    i2.  . . _ 

<€.>=<€.       >  +    <6e.    >  .  (45, 

J  J  J 

_2  _2 

Assuming     (fiy     )    =     (6y  )    and  using  equation  (44), 

J 


<6€l2>  =    <(£>*    TbY^    )  ~-   Eai     •  (T    <6y    ».  (46) 

J  *•       J  J  *•       J 


where  covariance  components  of  inaccuracy  are  neglected  for  simplicity. 
It  must  be  borne  in  mind  that  the  correlated  components  of  accuracy, 
which  are  in  general  unknown,  may  be  of  primary  importance,    particu- 
larly in  long  term.     Then  equation  (45)  will  be  written  as 


<0    =     <€°lZ>  +  Ea1;-    (T     <6y    »•  (47) 

J       J 


Some  correction  must  be  made  for  the  matrix  equations  derived 
in  section  2.      But  since  the  second  term  in  equation   (47)   gives  the  same 

effect  as  white  noise  FM,    no  correction  of   a    factors  for  white  noise  FM 

i 

is  necessary.      With  the  particular  cases  for  flicker  noise  FM  shown  in 
figure   10  where    a     factors  are  constant,    no  correction  is  necessary 
either. 

From  equation  (33),    the  ensemble  time  error  and  its  mean  square 
will  be  corrected  to 


\*  =  EEwlfl  =  Y,Ywl(*01  +  s^1)  (48) 

M        V^      l     i  4"^     i     j  i 

1    j       j    j  x     j      j     j  j 
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and 


<t6  =   <(EE»ifi)2>+  <(EEwi6€i)2>  (49) 


V  -   (iEE»0).+  (iEIw>()') 

ijJJ  ^jJJ 

respectively. 

In  order  to  get  a  rough  numerical  estimation  of  the  "accuracy, 
we  assumed  constant  weighting  factors  for 


w.  =     1/N  =  const,    and   a„    ~    a 
1  I.  I 

J 


used  to  obtain  equation  (50)  as  the  mean  square  error  of  the  ensemble 
time.     Then: 


<T^>    =    ^"<(E^)2>    +   M'A.r2(6y2>  (50) 

where 

2 
A    =  Va„    (  <  1)     =    const. 


The  first  term  in  equation  (50)  is  the  noise  dependent  one  and  the 
second  is  "accuracy"  dependent,    assuming  constant  mean  square  of  the 
inaccuracy  for  each  interval.     The  mean  square  accumulated  error  for 
one  clock  will  take  the  form  of  a  product  of  the  Allan  variance  and  some 
coefficient,    so  if 


<(Zcb2>  =   r2a2(r)  k-M  , 

J  J         Y 
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then  eq  (50)  will  be  simplified  as 


<t^>/T2cr2(T)  =    k-M/N  +  MA(6y2)/a2(T)  (5i; 

My  y 


where    k    is  constant. 

For  the  generation  of  the  most  accurate  and  uniform  time  scale, 
the  first  term  in  eq  (51)   should  be  negligibly  small  compared  to  the 
second.      The  condition  for  that  will  be  given  by 


k-M/N    «MA(6y2)/(j2(T)  (52) 


for  which 


<T^>  «  MAT2(6y2)  •  (53) 


Notice  that  the  ratio  of  mean  square  frequency  error  due  to 
"accuracy"  to  the  Allan  variance  due  to  noises  is   related  to  this  con- 
dition.     In  eq  (5  3)  the  mean  square  ensemble  time  error  due  to 
inaccuracy  increases  linearly  with  M ,    as  does  the  prediction  error  as 
mentioned  before. 

Table  I  gives  a  numerical  comparison  between  the  error  dependent 
on  prediction  and  error  dependent  on  "accuracy,  "   related  to  condition 
(52).     If  one  has  a  clock   10  times  as  good  in  stability  as  a  primary  stand- 
ard in  accuracy,    then  in  seven  days  the  value  of  the  right  hand  in  eq  (52) 
will  be  700.      There  is  also  shown  in  the  table  numerical  examples  of  the 
coefficient    k  •  M    by  approach  B.      Take   (2)  for  instance,    namely  two 
calibrations  in  a  week  and  24-hour     calibration  duration:    as  may  be  seen, 
those  values  of   k  •  M    are   sufficiently  small  compared  to  700.      So  in  this 
case,    we  can  construct  a  good  time   scale  with  only  one  clock,    which  is 
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very  simple  for  the  generation  of  a  time  scale.  With  the  value  coming 
from  a  condition  of  two  calibrations  per  week  and  six-hour  calibration 
duration,  one  needs  several  clocks  for  condition  (52),  even  using  such 
good  clocks. 

5.     COMPARISONS  BETWEEN  TIME  SCALES 

We  have  obtained  the  result  that  the  mean  square  of  the  predicted 
time  error  with  a  primary  standard  increases  with  M  (or   t    ,    elapsed 
time)  for  white  noise  FM,    flicker  noise  FM,    or  primary  standard 
inaccuracy  due  to  random  uncorrelated  noise.     One  can  see  this  from 
the  sample  shown  in  figure   18  by  comparing  curves  A  and  B  with  curve 
W    which  is  time  error  due  to  pure  white  noise,    which  also  represents 
time  error  due  to  "inaccuracy":    those  curves  appear  to  have   similar 
slopes . 

It  is  well-known  that  the  mean  square  of  a  simply  accumulated 

time  error   (without  prediction  and  calibrations)  increases  proportionally 

2 
to  elapsed  time    t      for  white  noise  FM  and  to    t      for  flicker  noise  FM. 

For  large  value  of  t      (^    10  days,    for  example),    however,    flicker  noise 

FM  will  be  usually  dominant  (see  fig.    17),    and  for  larger  value  of  t    , 

_2 

the  dominant  spectral  density    S    (f)    may  be  proportional  to    f         for  which 

Y  3 

mean  square  time  error  is  proportional  to    t    . 

Curve  S  in  figure   18  which  is  a  time  error  for  one  clock  can  be 


compared  for  that  with  curve  A  or  B  or  W.      (Notice  that  even  in  the 

-2 
prediction  method,    if  time  error  due  to  a  certain  noise,    f        for  exampl 

is  not  included  in  the  prediction,   then  it  will  become  dominant  and  give 

3 
the  same  time  error   (i.  e.  ,    t    )  as  its   simply  accumulated  time  error  af 

long  time  has  elapsed.  ) 
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Let  us  consider  the  ideal  case  where  the  predicted  time  error  is 
much  smaller  than  that  due  to  inaccuracy;    then  from  eq  (5  3),    the  mean 
square  time  error  is  given 

<T^>  ~   t.T<6y2)  (54) 

M    pred  1  } 


where   t      =  Mr.      The  mean  square  of  a  simply  accumulated  time  error 


averaged  over    N     clocks  will  be   [10] 


<T2.)    .  «   t2a2(t    )/N  (55) 

M    simple  1     y      1 


■where  constant  weighting  factors  over  clocks  are  assumed. 
Then  the  ratio  of  the  mean  square  errors  is 


(I2   )  J  (t2   )    .  =  r    N/M  (56) 

M    pred  M    simple  a 


where    r      =     <6y    )/<y     (t_), 
a  v      1 


-2,,     2 
r 

y 


2 

Let    r      =10      for  instance,    then  the  both  mean  square  errors  will 
a 

be  the  same  in  300  days   (^  1  year)  with   N  =  3    and   1000  days   (^  3  years) 
■with   N  =   10.      After  that  a  time   scale  generated  by  the  prediction  method 
-with  a  primary  standard  will  be  getting  better  and  better  relative  to  the 
other.      It  also  has  to  be  considered  that  the  "accuracy"  of  primary 
standards  will  continue  to  be  improved    [12]. 
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6.     ON  AN  APPROACH  USING  TIME  DIFFERENCE  DATA 
REFERRED  TO  THE  ENSEMBLE  TIME 

A  time  scale  generated  by  the  prediction  method  with  a  primary- 
standard  is  restricted  in  its  uniformity  by  the  accuracy  of  the  primary 
standard  used,    which  will  be  commonly  worse  than  the  frequency 
stability  of  clocks.     The  uniformity  of  a  time  scale  by  the  method  of  a 
simple  weighted  sum  over  clocks,    namely  without  prediction  and  cali- 
bration,   may  be  much  improved  for  a  short  time  range  by  increasing 
the  number  of  clocks.      (See  eq  (55).  )     So,    it  will  be  natural  to  expect  to 
obtain  a  better  uniformity  and  accuracy  in  a  time  scale  by  generating 
the  time  scale  with  the  optimum  time  prediction  and  with  an  initial 
calibration,    but  without  intermediate  calibrations. 

In  this  case,    one  might  use  the  ensemble  time  as  the  reference, 
instead  of  a  primary  standard.     So  the  time  difference  (or  frequency) 
data  of  each  clock  used  for  the  frequency  prediction  with  approach  A 
would  be  referred  to  the  ensemble  time  as   shown  with  notation   t   in 
figure   19.     The  following  discussion  on  this  approach,    however,    will 
show  that  little  improvement  can  be  expected. 

Let  the  predicted  frequency  or  time  difference  be  given  as  a 
linear  function  of  time  difference  data--L  in  number: 


T&  =  flM(AtM-r  ••■•  AtM-L> 


L 

i         .    i 

LMi         M  -  I 

1  =  1 


(57; 


T, 


where    r.    is  the  length  of  each  (prediction)  interval,    f ,  ,    is  a  linear 

M  B  u  ■  M 

operator  on    {At.}  at  the  M       interval  and    At  are  time  difference 
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data  for  i       clock  referred  to  the  ensemble  time  at  M  -  I        interval 

whose  duration  is    T  ,         (see  fig.    19);. 

M  -  I 


111  11 

At  =    t  -  t  and      t  =  x  -  T 

M-l  M-l        M  -  I  -  1  M-J&  M-l        M-J&-1 


so 

i  i  i  "  A 

M  -  I  =     (XM  -  I  '  XM  -  I  -  l]   '   (1M  -  I  '  TM  -  I-  1 


Ax1  -     AT  .  (5  8) 

M-  I  M  -  I  v 


In  the  approaches  so  far  mentioned,    the    a  =    a      in  eq  (57),    so 

1  V J.  JJ»  Xj 

that   f       =   f     (see  eq  (2)). 


The  predicted  time  error  may  be  expressed  as 


€M    =   XM  "  XM    =    AxM  "    TyM    =    AXM  "  fM(AxM  -  i  +    fM(ATM  -  i  ' 

(59) 

The  first  term  in  (59)  is  the   real  time  increase  during  the  M       interval 
and  the  second  is  a  prediction  term  for  the  first.      The  third  term  came 
from  using  ensemble  time  referred  data,    and  is,    in  a  sense,    redundant. 
Using  eq  (59),    the  ensemble  time  error  will  be  expressed  from  eq   (33), 
assuming  fixed  weighting  factors  for  each  interval  as 


.  M  ... 

t  ^yw'rtAx'-f1^1   j  +  f!(Ai.     n 

m     Z^     Z^        j      j       y-i         jv      j-r; 

i         j  =  l 


(60) 
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where    w    =   w     for  all   j    is  assumed.     Again  the  first  term  is  the  real 

J 
time  increase  during  the  0-M  interval  and  the  rest  is  a  prediction  for 

the  first. 

Now  let  us  reduce  the  third  term.     Similarly  to  eq  (60),    we  will 

get  the  following  two  equations  and  so  an  equation  for  their  difference 


L  y-Jj  i_         i.  L 

i  j=l         Jl       Jl        Jl 

1         Jl 


T.    ,       V  w1  V  {Ax.1  -  f.^Ax/     J+f^AT.       .)} 


and 


i     J-A-l  i  i  i  i 

\  J!=l  111  11 


SO 


^T._4    =   ^''(ix^  -  f^i-,/  +    i^i,-,/?  •  <61) 


Applying  the  linear  function    f     and  next     Sw  ,    we  will  obtain 

J  i 


11  A 


(62) 


+  F.-  F.       (AT.    ,     „   . 
J       3-1         J-A-^. 


where    F      =      T  w1-  f1    and  may  be  called  an  average  linear  function  over 

J  i  J 

clocks. 
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Substitution  of  eq  (62)  into  eq  (60)  gives 
M 


i     j=i 


)  +    F.(Ax'       ) 


(63) 
M 


^•f1      (Ax1  )l+y^F    •  F       (AT 

j      j  j  j  j=i     J        J  1 


which  again  includes  the  last  term- -a  function  of  ensemble  time 
difference. 

Repeating  the  above  procedure  and  defining  the  error  function 


Af1     =     f.1    -   F      for  all   k,  (64) 

K.  K.  ic 


we  will  finally  obtain  eq  (65), 
M 


Tw  =  V  w  V  TAx    -   Af  (Ax     J  -  F.-  Af     .(Ax1    n   n  ] 
i         j  =1  i 


F.-  F.       •   Af1    n    n    (Ax1    n    A      n   ) 


F  •  F        •  F  .  .  .  F  •  Af1  (Ax1  ) 

j    i-i    j-i-ij       j-i-ij.....^^     J-*-^..--V-3     J-£^r---V-2 

F  •  F        •  F  F  •  f 1  (Ax1  ) 

j    m    j-i-jj      j-i-i,...-^^  J-^-^----V_2     '-"r-ViJ 


M 
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where 


Ri  i  =  Fi'Fi-r"  Fi  it        i      (ATi  a        i      >-*°-    (66) 


The  last  two  terms  consist  of  real  time  difference  data  belonging 
to  zero  or  before  zero  intervals.     The  last  term   R    is  a  function  of 
initial  values  of  the  ensemble  time  error,    so  it  seems  reasonable  to 
assume  this  zero.     If  the  linear  functions  for  all  the  clocks  are  all  the 
same,    i.  e.  ,    if  the  clocks  are  statistically  identical,    then  all  the  error 
functions  in  eq  (64)  are  equal  to  zero,  then  all  the  medium  terms  in 
eq  (65)  will  disappear  to  give 


M 

(67) 

The  important  thing  to  notice  in  this  equation  is  that  the  ensemble 

*  i  i 

error    T  ,    does  not  depend  on  the  clock  values    (x    ,    .  .  .  ,  x,  ,    , )  actually 
M  ^  1  M-l  ' 

occurring  in  the  interval,    but  only  on  values  predicted  for  them  from 

values  occurring  before   j  =  1.      Thus,    no  new  information  is  added  in 

*  ens 
forming  the  estimate    T  based  on  clock  frequency  estimation  using 

a ens  a  ens 

the  ensemble  values    (T         ,    ...,T  ).     If,    in  fact,    the  clocks  started 

at  j  =  0,    (x       =0    for   m  ^  0),    eq  (67)  is  merely  the  simple  weighted 
m 

average  without  prediction 


T      =     >     w   (x , ,  -  x0 
M         jLj         v    M         0 


i  =  l 
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In  either  case,   the  results  of  using  approach  A  shown  in  figure  3  show 
that  for  flicker  noise  FM,    the  error    T        will  grow  rapidly,    at  least  as 
fast  as  \M/N,    as  we  advance   M   units  from  the  last  calibration. 


It  is  easily  shown  that  this  result  does  not  depend  on  the  fact  that 
we  have  weights    w      =  w  .      When  the  weights  are  allowed  to  vary,   the 
clocks  remaining  statistically  identical,    we  obtain  instead  eq  (67) 


N    M 

1=1  J  =  l 


F.-F.    ..  .  .  F.    „  „  f.    . 


Ax1 


M-2 


j-i 


)} 


M-l 


where  again  the   second  (correction)  term  does  not  depend  on  clock 
values  which  occur  for    j   >  0,    and  vanishes  for  clocks   started  at   j   =  0, 

In  the  general  case  where  the  clocks  are  required  to  be 
statistically  identical,    we  obtain  by  the  same  procedure 


N     M 


M  £- i   ^— I   I  11  &  1  1  - . 

i=l  j=l       J  J         J       J 


F.-Ag1    .(Ax1-! 
J        J-i      J 


F.-F       •  Ag1  (Ax1         „  ) 

J       j-i         ^j-i-ij  j-i-Y^^ 


F.-F.       .  .  .    F 
J       J-i  J-i 


-F.-F. 
j       J-i 


M 


R. 


j  =  l 


j,i 


^gi 


M-4 


M-3 


j-i 


i 
'j-i 


(A: 


M-3 


>i 


I 


(Ax: 


M-2 


j-i 


-i 


) 


M-2 


-L 


M 


i»] 


(68) 
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where 


(69: 


N 
r         =     7       w     I       .  Sf         =     w     f 

k        La      k  k  '  gk  k  k 

i  =  I 


N 

A%  =  wk  (£k  -  E  wkfk» 

s  =  l 


Equation  (68)  does  contain  prediction  terms  dependent  on  clock 

values  occurring  on  the  interval   j   =   (1,  M  -  1),    so  that  new  information 

is  being  added.      But  it  does  not  seem  likely  that  the  error    T        would 

M 

show  significantly  reduced  growth  rate  from  that  obtained  in  the  case 
with  statistically  identical  clocks.      It  is  important  to  note  that  eq  (68) 
for    T        is  nonlinear   (of  degree  M)  in  the  weights    w      ,     so  that  the 
optimum  estimation  of  these  weights  is  very  difficult. 

The  conclusion  appears  to  be  that  while  frequency  calibrations 
against  a  primary  standard  can  greatly  improve  the  accuracy  and 
stability  of  the  time  scale,    for  the  noise  processes  considered  in  this 
paper,    no  improvement  in  these  quantities  can  be  gained  (at  least  in 
the  case  of  statistically  identical  clocks)  by  using  the  time   scale  itself 
as  a  primary  standard. 
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7.     CONCLUSIONS 

We  have  obtained  the  result  that  one  can  construct  with  a  small 
number  of  clocks  a  time  scale  whose  accuracy  is  limited  by  the  accu- 
racy of  a  primary  standard.      We  have  also  shown  that  the  time 

dispersion  of  the  time  scale  thus  calculated  may  be  proportional  to  the 

2 
square  root  of  the  elapsed  time    t    in  a  statistical  sense  instead  of   t 

for  the  time  scale  without  intermediate  calibrations,    so  the  former 
may  become  much  better  in  both  accuracy  and  uniformity  than  the  latter 
after  a  long  time  has  elapsed.      Furthermore,    algorithms  using  a 
primary  standard  have  the  advantage  that  it  is  easier  to  identify  fre- 
quency drift,    jumps,    or  aging  of  the  clocks  used. 

It  must  be  emphasized  that  these  conclusions  are  based  on  clocks 
having  only  continuous  random    f      and   f         FM  noise  with  constant 
amplitudes.      In  practice,    these  amplitudes  can  vary  in  time,    and  many 
other  anomalous  forms  of  behavior  occur,    which  will  strongly  influence 
the  choice  of  both  the  clock  ensemble  and  the  time  scale  algorithm  to  be 
used. 

We  have  shown  that  (at  least  for  statistically  identical  clocks), 
no  advantage  is  gained  by  using  the  time   scale  itself  as  a  "primary- 
standard"  for  intermediate  calibrations. 

It  is  also  true,    however,    that  the  time  scale  averaged  over  a 
large  number  of  clocks  without  intermediate  calibrations  may  have 
a  much  better  uniformity  for  a  short  time  range  than  that  depending 
on  a  primary  standard  utilized  as  in  the  above  methods  of  prediction, 
where  the  short  term  uniformity  of  the  time  scale  is  limited  by 
accuracy  of  the  primary  standard.      Therefore,    it  may  be  possible  to 
generate  the  two  kinds  of  time  scale  simultaneously  and  actually 
utilize  the  former  time  scale  periodically  calibrated  by  the  latter. 
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Fig.  1  (a)  Linear  drift  of  time.  (b)  Time  error,  x  ,  due  to  random 
noises  and  its  prediction  error,  e  .  (c)  A  group  of  predicted 
time  errors,      e1.    and  the  ensemble  time,    T. 
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Fig.    2     Approach  A. 
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Fig.    3       Mean  square  prediction  error  by  approach  A  as  a  function  of 
position  of  a  calibration,     f      and   f~*    noise  FM  refer  to  white 
and  flicker  noise  FM,    respectively,    and   f    +   f~      noise  FM 
refers  to  a  mixed  noise  of  both  with  power  ratio    r  =  1. 


42 


2.0 


h  =  i 


K> 


CM 


<\J 


I-    1.5 


E 

CM 


1.0 


L=2 


f 


L  =  5 


r 


(ex)^s=2 


M  I   M  I   MTI 

CAL      CAL     CAL 

(FIXED) 


/PREDICTION 
\       INTERVAL 


0     12     3     4     5     6 
Spacing  of  Calibration 


Fig.  4  Mean  square  prediction  error  by  approach  A  as  a  function  of 
spacing  of  calibration.  fu  and  f "  *  refer  to  white  and  flicker 
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refer  to  white  and  flicker  noise   FM,    respectively. 
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Fig.    7       Mean  square  prediction  error  by  approach  B  as  a  function  of 
position  of  a  calibration.      Calibration  duration  is  one  day. 
f      and   f "  *■   refer  to  white  and  flicker  noise  FM,    respectively. 
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Mean  square  prediction  error  by  approach  B  as  a  function  of 
number  of  calibrations  in  a  week.      Calibration  duration  is  one 
day.     f      and    f"       refer  to  white  and  flicker  noise  FM,    respec- 
tively.     The  dashed  continuation  for  the    f         curve   shows  the 
case  where  the  optimum  positions  of  calibrations  are  selected 
for    L   =  4  ~   6;  Saturday  and  Sunday  are  used  for  calibration 
days. 
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Mean  square  prediction  error  by  approach  B  as  a  function  of 
calibration  duration,      f      and    f~       refer  to  white  and  flicker 
noise  FM,    respectively. 
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Fig.    10        Optimum  filter   response  functions  for  flicker  noise  FM  by- 
approach  B.      (a)    and  (b):     Two  calibrations  per  week, 
(c):    Three  calibrations  per  week,      (d):    Five  calibrations  per 
week.     Calibration  duration  is  one  day. 
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Fig.     11        Mean  square  prediction  error  by  approach  B  as  a  function  of 
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refer  to  white  and  flicker  noise  FM,    respectively. 


50 


CM 


h=l  L=l  a  2Cals/Week 


i-     4 


2- 


u- 

0 


J I I L 


14 


21 


28 


M 


Fig.    12        Mean  square  of  accumulated  prediction  error  by  approach  A 
(L  =   1).     Calibration  duration  is  one  day.     f     and   f 
to  white  and  flicker  noise  FM,    respectively. 


refer 


51 


h=l   L=2  a  2Cals/Week 


OJ 


CM 


X/ 


Fig.    13        Mean  square  of  accumulated  prediction  error  by  approach  A 
(L  =  2) .     Calibration  duration  is  one  day.     f      and   f~       refer 
to  white  and  flicker  noise  FM,    respectively. 
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Fig.    14        Mean  square  of  accumulated  prediction  error  by  approach  A 

as  a  function  of  number  of  calibrations  in  a  week.     Calibration 
duration  is  one  day.     f  ,  f        and   f  +  f        refer  to  white  noise 
FM,    flicker  noise  FM,    and  a  mixed  noise  of  both,    respectively. 
The  dashed  continuations  of  the  curves  show  the  case  where  the 
optimum  positions  of  calibrations  are  selected  for  the  horizontal 
axis    4  ~  7;  Saturday  and  Sunday  are  used  for  calibration  days. 
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white  noise  FM,    flicker  noise  FM,    and  mixed  noise  of  both, 
respectively. 
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A  numerical  comparison  between  predicted  time  error  by  approach  B 
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ly). Available  only  from  NBS  Boulder  Labora- 
tories. Ordering  and  cost  information  may  be 
|  obtained  from  the  Program  Information  Office, 
National  Bureau  of  Standards,  Boulder  Colo- 
rado 80302. 


NONPERIODICALS 

Applied  Mathematics  Series.  Mathematical  tables, 
manuals,  and  studies. 

Building  Science  Series.  Research  results,  test 
methods,  and  performance  criteria  of  building  ma- 
terials,  components,   systems,   and  structures. 

Handbooks.  Recommended  codes  of  engineering 
and  industrial  practice  (including  safety  codes)  de- 
veloped in  cooperation  with  interested  industries, 
professional  organizations,  and  regulatory  bodies. 

Special  Publications.  Proceedings  of  NBS  confer- 
ences, bibliographies,  annual  reports,  wall  charts, 
pamphlets,  etc. 

Monographs.  Major  contributions  to  the  technical 
literature  on  various  subjects  related  to  the  Bureau's 
scientific  and  technical  activities. 

National    Standard     Reference    Data    Series. 

NSRDS  provides  quantitative  data  on  the  physical 
and  chemical  properties  of  materials,  compiled  from 
the  world's  literature  and  critically  evaluated. 

Product  Standards.  Provide  requirements  for  sizes, 
types,  quality,  and  methods  for  testing  various  indus- 
trial products.  These  standards  are  developed  co- 
operatively with  interested  Government  and  industry 
groups  and  provide  the  basis  for  common  understand- 
ing of  product  characteristics  for  both  buyers  and 
sellers.  Their  use  is  voluntary. 

Technical  Notes.  This  series  consists  of  communi- 
cations and  reports  (covering  both  other-agency  and 
NBS-sponsored  work)  of  limited  or  transitory  interest. 

Federal  Information  Processing  Standards 
Publications.  This  series  is  the  official  publication 
within  the  Federal  Government  for  information  on 
standards  adopted  and  promulgated  under  the  Public 
Law  89-306,  and  Bureau  of  the  Budget  Circular  A-86 
entitled,  Standardization  of  Data  Elements  and  Codes 
in  Data  Systems. 

Consumer  Information  Series.  Practical  informa- 
tion, based  on  NBS  research  and  experience,  cover- 
ing areas  of  interest  to  the  consumer.  Easily  under- 
standable language  and  illustrations  provide  useful 
background  knowledge  for  shopping  in  today's  tech- 
nological marketplace. 

CATALOGS  OF  NBS  PUBLICATIONS 

NBS  Special  Publication  305,  Publications  of 
the  NBS.  1966-1967.  When  ordering,  include 
Catalog  No.  CI 3.1 0:305.  Price  $2.00;  50  cents  addi- 
tional for  foreign  mailing. 

NBS  Special  Publication  305,  Supplement  1, 
Publications  of  the  NBS,  1968-1969.  When  order- 
ing, include  Catalog  No.  C13.10:  305/SuppI.  1.  Price 
$4.50;  $1.25  additional  for  foreign  mailing. 

NBS  Special  Publication  305,  Supplement  2, 
Publications  of  the  NBS,  1970.  When  order- 
ing, include  Catalog  No.  C13.10:  305/Suppl.  2.  Price 
$3.25 ;  85  cents  additional  for  foreign  mailing. 


Order    NBS    publications    (except     Bibliographic    Subscription    Services) 
horn:   SuP£nntendent  of  Documents.  Government  Printing  Office.  Wash- 
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